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ABSTRACT 

The migration of a planet through a gaseous disc causes the locations of their resonant 
interactions to drift and can alter the torques exerted between the planet and the 
disc. We analyse the time-dependent dynamics of a non-coorbital corotation resonance 
under these circumstances. The ratio of the resonant torque in a steady state to the 
value given by Goldreich & Tremaine (1979) depends essentially on two dimensionless 
quantities: a dimensionless turbulent diffusion time-scale and a dimensionless radial 
drift speed. The dimensionless diffusion time-scale is a characteristic ratio of the time- 
scale of turbulent viscous diffusion across the librating region of the resonance to the 
time-scale of libration; in the absence of migration, this parameter alone determines 
the degree of saturation of the resonance. The dimensionless radial drift speed is the 
characteristic ratio of the drift speed of the resonance to the radial velocity in the 
librating region; this parameter determines the shape of the streamlines. When the 
drift speed is comparable to the libration speed and the viscosity is small, the torque 
can become much larger than the unsaturated value in the absence of migration, but 
is still proportional to the large-scale vortensity gradient in the disc. Fluid that is 
trapped in the resonance and drifts with it acquires a vortensity anomaly relative 
to its surroundings. If the anomaly is limited by viscous diffusion in a steady state, 
the resulting torque is inversely proportional to the viscosity, although a long time 
may be required to achieve this state. A further, viscosity-independent, contribution 
to the torque comes from fluid that streams through the resonant region. In other 
cases, torque oscillations occur before the steady value is achieved. We discuss the 
significance of these results for the evolution of eccentricity in protoplanetary systems. 
We also describe the possible application of these findings to the coorbital region and 
the concept of runaway (or type III) migration. 

Key words: accretion, accretion discs — galaxies: kinematics and dynamics - 
hydrodynamics — planets and satellites: general 



1 INTRODUCTION 

Corotation resonances play an essential role in the grav- 
itational interaction between a planet and the disc in 
which it forms. The perturbing potential associated with 
a planet (or other satellite) orbiting in a disc can be de- 
composed into a series of Fourier components with differ- 
ent azimuthal wavenumbers m and angular frequencies ui 
JColdreich fc Tremaind ll980H . A fluid element in the disc 
at radius r, where the angular velocity is fi(r), experi- 
ences a perturbing force with Doppler-shifted frequency ui = 
uj — rafi. Resonant responses to horizontal forcing occur at 
radial locations where ui = (the corotation resonance) and 
£b = ±k (Lindblad resonances), where n(r) is the epicyclic 



frequency in the disc. Under certain assumptions, includ- 
ing a linear approximation, simple formulae can be obtained 
for the localized resonant tor ques exerted at these locations 
jGoldreich fc Tremainelll979l hereafter GT79). The associ- 
ated exchanges of angular momentum and energy cause the 
orbital parameters of the planet, in particular its semimajor 
axis and eccentricity, to evolve, and also play a role in shap- 
ing the mass distribution in the disc. These processes are 
of fundamental importance in determining the properties of 
planetary systems. 

In a linear analysis the torque exerted at a corotation 
resonance is found to be proportional to the local radial gra- 
dient of vortensity (i.e. vertical vorticity divided by surface 
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density). However, the perturbed flow in the corotation re- 
gion has islands of libration within which fluid is trapped and 
corotates, on average, with the potential. The mixing within 
these islands tends to erase the gra dient of vortensity and 
reduce the corotat i on torque to zero iGoldreich fc Tremaind 
ll98ll:lLubowlll990l : IWardll99ll.ll992ri . The resonance is then 
said to be completely saturated. Since viscous diffusion tends 
to reestablish the vortensity gradient, the level of saturation 
achieved in a steady state depends on the strength of the 
forcing and the viscosity. This behaviour is in contrast to 
that of Lindbla d resonances, where th e torque is insensitive 
to nonlinearity jYuan fc Cassen|ll994h and can be reduced 
significantly only by removing mass from the resonant re- 
gion. 

In the case of a planet with a circular orbit, the angular 
pattern speed Sl p = w/m of each non-axisymmetric poten- 
tial component is equal to the angular velocity of the planet, 
and the corotation resonances are all coorbital, occurring at 
(very nearly) the same radial location as the planet's or- 
bit. The coorbital corotation resonance is difficult to treat 
analytically because of the simultaneous presence of reso- 
nant components of all azimuthal wavenumbers and the high 
degree of nonlinearity. It is also in this region that three- 
dimensional effects are most important. The presence of a 
gap or partial clearing around the planet's orbit, together 
with shocks and accretion stream s, complicates the problem 
furthe r. An important study by iBalmfort h fc Korvcanskvl 
(2001) presented an asymptotic reduction of the coorbital 
region, allowing for a modest degree of nonlinearity, to a 
more tractable problem requiring the numerical solution of 
a partial differential equation. The time-evolution of the so- 
lution shows the formation of transient vortices followed by 
a mixing of vortensity across the corotation region, resulting 
in a saturation of the corotation torque. 

Some recent direct numerical simulations of the in- 
teraction between a mobile planet and a disc have drawn 
attention to the importa nce of the coorbital region. 
iMasset fc Papaloizoul (|200.?t) found that very fast migration 
of a Saturn-mass protoplanet could occur in a disc somewhat 
more massive than the minimum-mass solar nebula. They 
presented a simplified analytical model allowing them to es- 
timate the torque exerted by gas that crosses the planet's 
orbit as the planet migrates radially through the disc. Since 
this torque is proportional to the migration rate, it mod- 
ifies (and typically reduces) the effective inertial mass of 
the planet, allowing it to migr ate faster in response to non - 
coorbital torques. In addition. IMasset fc Papaloizoul (120031) 
suggested that, if the effective inertial mass became negative 
and a time delay also occurred in its determination, the ra- 
dial motion could be unstable and give r ise to very fa s t (and 
generally inward) 'runaway' migration. lArtvmowiczl <l2004l) 
presented a preliminary model for the torque due to the 
asymmetry in the librating region owing to planet migra- 
tion, based on the behaviour of particle orbits. These ef- 
fects have been described as 'type-Ill' migration to differ- 
entiate this regime from the previously identified regimes 
of migration of e mbedded plan ets (type I) and gap-opening 
planets (type II) <Wardlll99'fl . Rece ntly the numerical re- 
su lts of [Masset fc Papaloizoul l|2003f) have been challenged 
bv ID ' AngeTo"F3ate"fcHibowr(l2005l) . who found that the ef- 
fect is strongly or completely suppressed when the region 



close to the planet is better resolved. The reality of runaway 
migration is therefore in some doubt. 

In the case of a planet with an eccentric orbit, poten- 
tial components are present at first order in the eccentricity 
that have angular pattern speeds different from the mean 
motion of the planet. The principal effect of the associated 
corotation and Lindblad torques is to cause the eccentric- 
ity to evolve. Provided that a deep gap is cleared around 
the planet, so that the coorbital eccentric Lindblad reso- 
nances are ineffective, there is a fine balance between the 
growth of eccentricity through eccentric Lindblad resonances 
and its decay thro ugh eccentric corotation resonances, whic h 
are non-coorbital iGoldreich fc Tremaineil980tlWarol l988L 
if many such resonances are able to compete with each 
other. Even a partial saturation of the corotation torques 
therefore promotes the growth of the planet's eccentricity. 
It should be noted, however, that whether the eccentricity 
grows or decays depend s on a number of additional factors. 
Gold reich fc Saril J2003ft discussed some of the relevant is- 
sues but further work is required to treat the development 
of eccentricity in the disc and the conservative secular ex- 
change of eccentricity between the planet and the disc. 

We recently analysed the saturation o f the non-coorbital 
corot ation resonance in a gaseous disc jpgilvie fc Lubow 
|2003l here after Paper I). Unlike the ca lculations "of lMassel 
j200l[) and lMasset fc Papaloizoul (l2003l> . which, owing to the 
extreme complexity of the coorbital region, necessarily have 
the nature of a toy model, our calculation is based on an 
asymptotically exact reduction of the problem making min- 
imal assumptions. We showed that the steady corotation 
torque for any non-coorbital resonance is reduced below the 
value specified by GT79 by a factor t c (p) depending on a sin- 
gle dimensionless parameter p, which measures the strength 
of the forcing relative to the effects of viscosity in the disc. 
We computed the function t c {p) numerically and derived 
analytical a pproximations for s mall and large p. Our result 
was used bv lGoldreich fc Saril i2003T) in their theory of the 
growth of the eccentricity of a protoplanet through a finite- 
amplitude insta bility, and ou r findings have recently been 
corroborated by M asset fc Og ilvie ( 2004) using localized nu- 
merical simulations. 

The principal aim of this paper is to extend the anal- 
ysis of Paper I to include the effect of planetary migration. 
We also examine the time-dependent approach to a steady 
state and emphasize the interpretation of the solutions in 
real space rather than Fourier space. Finally, we attempt 
to relate our findings in a preliminary way to the coorbital 
region and the issue of runaway migration. 



2 INCLUSION OF TIME-DEPENDENCE 

In Paper I we carried out a systematic asymptotic analysis of 
the corotation region of a three-dimensional, barotropic, vis- 
cous, non-self-gravitating disc subject to a uniformly rotat- 
ing external potential perturbation. To calculate the steady 
torque exerted on the disc, we sought a solution that is 
steady in the corotating frame of reference. We made use 
of the small parameter e, which is a characteristic value of 
the angular semithickness H/r of the disc. In units such that 
the corotation radius and the corresponding angular veloc- 
ity are of order unity, we adopted natural scalings such that 
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the width and height of the corotation region are O(e). We 
introduced scaled radial and vertical coordinates, t and £, 
to resolve the inner structure of this region. 

It is a simple matter to restore time-dependence to the 
original problem. The characteristic time-scale for establish- 
ing the steady solution is 0(e _1 ), this being typical of both 
the libration time-scale and the viscous diffusion time-scale 
across the corotation region under our scaling assumptions. 
We therefore allow the solution to depend on a 'slow' time 
variable r = et, and the effect is to replace each instance of 
the operator 



as in equation (17) of Paper I, with 



3_ 

dr 



d_ 
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A further adaptation allows for a 'very slow' time- 
dependence of the external potential. This is natural if we 
consider a Fourier component of the tidal potential of a 
planet that migrates through the disc on a time-scale that is 
long compared to the orbital time-scale and also compared 
to the characteristic time-scales of libration and viscous dif- 
fusion across the corotation region. In this case both the 
radial structure and the angular pattern speed of the poten- 
tial evolve very slowly in time. We allow for this formally by 
introducing a very slow time variable T — e 2 t and writing 
the potential perturbation in the midplane z = as 



where 



¥ = 



fip(T) AT, 



(3) 



(4) 



4> being the usual azimuthal angle in an inertial frame of ref- 
erence. The pattern speed of the potential is then S7 P (T), a 
very slowly varying angular frequency of order unity, and 
the corotation radius is r c (T), defined by the condition 
0(r c (T)) = fl p (T). 

We transform from spatial coordinates (r, <f), z) to 
C)i with 



t = 



Tc{T) 



c = -. 

e 



Using the chain rule we find 
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(6) 
(7) 
(8) 
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We proceed to solve the fluid dynamical equations as in 
Paper I by expanding the solution in powers of e. Whereas 
previously we looked for a strictly steady solution in the 
corotating frame, now the solution depends on the slow and 
very slow time variables r and T. The time-derivative asso- 
ciated with the dependence on T is small, 0(e 2 ), and does 
not affect the equations to the order that we considered. On 



the other hand, there is a new 'advective' time-derivative 
— (Ar c / At) (d/d£) at O(e) that does affect the analysis, along 
with the derivative with respect to r. The net effect is that 
the operator Q is replaced with 



d f , Ar c \ d d 
^ + { U °-AT)di +1 % 



(10) 



By analyzing the problem in a frame of reference that moves 
with the resonance (hereafter referred to as the comoving 
frame), we transform the planetary migration into a uni- 
form radial drift of the gas through the resonant region (cf. 
iMasset fc PapaloizodEooli) . 

When we remove the e-scalings and present the reduced 
equation for the enthalpy perturbation h'(x,<p,t) in dimen- 
sional form, we obtain 



d 
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(11) 



where we recall that x — r — r c is the radial distance from 
corotation, 



1 d&_ 
2r~B~dlp~ 



(12) 



is the leading-order radial velocity perturbation induced by 
the tidal potential 2B = (l/r)(d/dr)(r 2 fi) is the vertical 
vorticity of the unperturbed disc, c is a certain vertical aver- 
age of the sound speed, v is the mean kinematic viscosity, E 
is the surface density and D M = 91n(i/E)/91n E. In equation 
lilt all coefficients are to be regarded as independent of x 
and evaluated at r = r c , except in the one place where x 
appears explicitly. The quantities and v! are regarded as 
functions of ip only. The time-dependence of the coefficients 
is neglected, on the basis that the time-scale for establishing 
the solution is short (measured by r) compared to that on 
which the coefficients vary significantly (measured by T). 

As in Paper I, we consider a single potential component 
of the form 



$ = ^cos{rmp), 



(13) 



and the solution for h! may be assumed to have the same 
periodicity in ip. The total tidal torque exerted in the coro- 
tation region can be written as 



T c (t) 



rE 
"75" 



h'(x,<p, Ax Atp. 



(14) 



The boundary condition as \x\ — > oo is that dh' /dx — > 0. 
This ensures that the disturbance is localized in the corota- 
tion region in the sense that the vortensity gradient at large 
x is just that of the unperturbed disc. Generally, h 1 does not 
tend to zero as \x\ — + oo (although it is bounded) because 
the corotation torque must be balanced in a steady state 
by an adjustment of the viscous stress (and therefore the 
surface density) across the corotation region. 
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3 REDUCTION TO A DIMENSIONLESS 
FORM 

As in Paper I, we rewrite equation pip in a dimensionless 
form by means of the transformations 



dr\ 



* mc ( d< I ' 



c 

X -, 

K 



We then obtain 
d 



V = — , 
m 



dx 2 



(15) 
(16) 



where 
a = 2 



dx 2 

/ dlnr \ * 
\ dlnfi/ 
dft D u 
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u _-2rn^^ 
dr k 

1 / dr \ k 
^ml-dfij 

1 / dr\ k 2 dr c 
= m V~ dT2/ "?"dT' 
The ratio 
_ d _ 2rff dr c 
~ a ~ m<3> dt 



(17) 

(18) 
(19) 
(20) 
(21) 

(22) 



measures the relative effect of the new advective term asso- 
ciated with the drift of the corotation resonance. 

The corotation torque can be written in the form 



where t c is dimensionless, and 



-EX 



(23) 



(24) 



2(dfi/dr) dr \BJ 
is the torque formula of GT79. 

4 ANALYSIS IN REAL SPACE 

4.1 Interpretation as a vortensity equation 

In a two-dimensional, barotropic, inviscid flow, the vorticity 
equation takes the form 

D / uj. 



Dt V E ) 



= 0, 



(25) 



where lo z = e z ■ (V x u) is the (vertical component of) vortic- 
ity. The quantity lj z /E is sometimes referred to as the poten- 
tial vorticity or (among accretion disc theorists) the vorten- 
sity. It is not generally conserved in a three-dimensional flow, 
nor is it conserved in the presence of viscosity. 

We consider the related quantity Q — \tx(S/lv z ). In the 
unperturbed disc this is equal to Q = ln(S/S) — In 2 and 
has a purely radial gradient AQ/Ar — (d/dr) ln(E/B). This 
gradient appears in equation 11 I li and can be regarded as 
constant across the corotation region. In the perturbed disc, 
under our scaling assumptions, 



Q = Q + Q' = Q+% 



1 d 2 h' 
k 2 dx 2 



(26) 



Equation 111! cannot generally be written as a closed equa- 
tion for Q because of the viscous terms. However, a simpli- 
fication occurs when the dimensionless parameter 

6=-2rn^% (27) 
dr k 

equal to 3D M in a Keplerian disc, is equal to unity. As in 
Paper I, we adopt this convenient assumption, anticipat- 
ing that our results will not depend sensitively on it. When 
6=1, equation I I U can be divided through by k 2 and in- 
terpreted in the form 



« + ..v*. 

where 



9 2 



dx 2 



dr ' 



( , dr c \ dQ 



(28) 



(29) 



is the velocity field in the comoving frame correct to a certain 
level of approximation. 1 Equation l|28p is equivalent to 



(30) 



provided that the diffusive V 2 operator is considered to act 
only on the more rapid variation of the perturbation Q' in 
the radial direction. It should be noted that the quantity 
DQ/Dt does not contain any contribution from the drift; it 
is a Galilean-invariant quantity and is most easily evaluated 
in the non-moving frame. Evidently, equation 1301 derives 
from a version of potential vorticity conservation. Three- 
dimensional effects are absent for the flows under considera- 
tion, which are quasi-two-dimensional because the motion is 
predominantly horizontal and maintains a quasi-hydrostatic 
balance in the vertical direction. Only when b — 1 do the 
viscous terms have a simple closed form in terms of Q. 



4.2 Streamlines of the dominant motion 

It is of interest to plot the streamlines of the velocity field 
1251 . which are contour lines of the streamfunction 



X(x, tp) = cos(m^) - r^-ip ■ 



rdO. 2 
2dr~ X ■ 



(31) 



These are not the exact streamlines of the fluid but represent 
the dominant motion u — Vx x e z in the comoving frame. 
In the rescaled variables of Section |3J the streamfunction is 
proportional to 



d9+ ii 2 , 



(32) 



and the shape of the streamlines depends only on the di- 
mensionless parameter 

d 2rB dr c , . 

v = - = — — — — , (33) 

introduced previously, which measures the drift speed rel- 
ative to the characteristic libration speed. The streamlines 

1 Since the solution varies more rapidly in the radial direction 
than in the azimuthal direction, u r is in fact correct to 0(e 2 ) 
while Utp is correct to 0(e). 
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are plotted in Fig.^for the cases v = 0, 0.5, 1 and 2 (those 
for negative values of v can be obtained by a reflection in the 
i-axis). Stagnation points occur where x = and sin 9 = v, 
which yields two solutions in < 9 < 2tt if |d| < 1 and none 
if \v\ > 1. In the absence of drift, the streamlines circulate 
both interior and exterior to the resonant radius, but form 
symmetric librating islands centred on the resonance. The 
drift breaks the symmetry of the flow with respect to the x- 
axis. If the drift is not too fast (\v\ < 1) the librating islands 
remain but they become asymmetrical and are diminished 
in size; this opening allows gas to flow through from the 
outer disc to the inner disc (if v > 0). A fast drift (|«j > 1) 
destroys the librating islands and allows a free streaming 
across the corotation region. For a given strength of poten- 
tial there is therefore a critical migration rate that changes 
the topology of the flow and inhibits the libration normally 
associated with the corotation resonance. 



4.3 Behaviour of the vortensity perturbation 

Now equation 12811 can be interpreted as an advection- 
diffusion equation with a steady source term. The quantity 
Q' that is being advected and diffused is (minus) the frac- 
tional vortensity perturbation, and the source term derives 
from the velocity perturbation u' induced by the forcing. 
In the presence of viscosity a steady solution will be ap- 
proached in which the advective and diffusive terms balance 
the source term, i.e. 



u-VQ' 



d 2 



,dQ 
dr 



(34) 



In Appendix [XJ we analyse the solution of an equation of 
this type in the limit of small viscosity. The conclusion is 
that the perturbation builds up on the closed streamlines 
of the librating region to a value, proportional to v~ x , at 
which it is equilibrated by outward diffusion to the region of 
open streamlines. We apply equation l|A140 of Appendix 1X1 
by setting the source term equal to — u- VQ — (dr c /dt)dQ/dr 
and noting that the first contribution integrates to zero. Also 
V 2 x is replaced by the constant value d 2 x/dx 2 because the 
radial derivatives are asymptotically dominant. Therefore 



vr \ 



dr \ dr c dQ 
dSi) ~dt dr 



(x - Xs) 



(35) 



in the librating region in a steady state, where Xe is the 
streamfunction of the separatrix. In the region of open 
streamlines outside, 



Q' 



-u — I dA, 

dr 



(36) 



where the integral is taken along the streamline from far 
upstream, and A is the travel time along the streamline. 



4.4 Estimation of the torque in the low-viscosity 
limit 

The corotation torque 11411 can also be expressed in terms 
of Q' , noting that the last term in equation 1261 integrates 
to zero: 

f2n poo 



-rE 



yf**L 

dip 



dx dip 



2rBE 
dQ/dr 



,dQ 
dr 



dA. 



(37) 



In the low- viscosity limit, the torque in a steady state is then 
(again using the results of Appendix [XJ 



2rBE /dr, 



m 2 vrSdQ 
2rB dr 



dQ 
dr 

1 
2 



sin 9 dA 



dx, 



(38) 



where the first integral is over the librating region and the 
second is over the region of open streamlines. 

The separatrix bounding the librating region is the 
streamline passing through the saddle point at x = 0, 
9 — 9 S = 7r — arcsin v. The equation of the separatrix can be 
written in dimensionless form as 



[-( 


' dr V 


_rB » 


v~d«/. 



1/2 



X(6), 



where X(9) is defined by 
X 2 = cos 9 — cos 9 a + v(9 ■ 



(39) 



(40) 



and is positive in the interval 9\ < 8 < 9 S , where 9i is the 
other value of 9 at which the separatrix crosses the a;-axis. 
After some rearrangement, and remembering that there are 
m copies of the librating region, we find the ratio of the 
torque to the GT79 value in the low-viscosity limit to be 



2^/2 d 2 
it 2 a 1 / 2 !/ 

1 



+ - 



7T 



°l[X(9)fd9 



sin 9 dA 



dx, 



(41) 



where the second integral is over a single copy of the stream- 
ing region, and the dimensionless streamfunction and travel 
time are defined by 



X = X 



rc 2 ( dfi\ 
k 2 \ dr J 



A = A 



K f dr \ 
mc V dQ J 



(42) 



4.5 Interpretation of the torques 

The physical meaning of this analysis is that, if the disc has 
a large-scale vortensity gradient, fluid that is trapped in the 
corotation resonance and forced to drift with it attempts 
to preserve its original vortensity. A steady state is achieved 
when the vortensity anomaly is equilibrated by viscous diffu- 
sion to the fluid outside the librating region. If the viscosity 
is small, the torque comes mainly from the librating region 
(the first term in equation 1381 or equation I41|l but there is a 
further contribution, differing from the GT79 torque only by 
a factor of order unity, from the fluid that streams through 
the resonance (the second term in the same equations). 

The first integral in equation l|41|l is a dimensionless 
quantity of order unity that decreases monotonically from 
32^/2/9 w 5.028 to as 1 1? | increases from to 1 and the li- 
brating region shrinks (cf. Fig.0. Therefore the first contri- 
bution to t c scales as d for small v. In dimensional 
terms this contribution to the torque scales as 



3/2 



dlnfi 
dlnr 



Q, 2 i 



/drc\ 2 d_ / E\ 
\dt) dr\BJ 



(43) 
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v = ; 0.5 




Figure 1. Streamlines of the dominant motion in the comoving frame, for various values of the drift-to-libration ratio v = d/a. The 
separatrices are marked as bold lines. A librating region of closed streamlines exists for \v\ < 1. The variable 6 is shown over two periods 
for greater clarity; the actual number of periodic copies is m. The increment of streamfunction between neighbouring contours is twice 
as large in the lower panels as in the upper panels, implying that the velocities are typically larger. 



A simple order-of-magnitude argument for this result can be 
given as follows. In the absence of turbulent diffusion, the 
vortensity in the librating region remains constant as the 
planet migrates. However, diffusion limits the contrast in Q 
between the librating region and the background to AQ ~ 
— (dr c /dt)t I/ (d<5/dr), where t v = S 2 /u is the viscous time- 
scale across the libration region of radial width 5 for which 
5 2 ~ $/[rB(— dfi/dr)]. Then the torque in equation E51 
can be understood as ~ — (AQ)(Y,rS)(rB)(dr c /dt). This is 
only a fraction of order —AQ of the rate of change of orbital 
angular momentum of the trapped region. The reason that 
only a fraction is required is that the region occupied by the 



moving trapped fluid is replenished by fluid of the ambient 
vortensity. 2 

The angular momentum taken up by the vortensity 
anomaly is somewhat analogous to the kinetic energy de- 
veloped by the entropy anomaly in stellar convection. Both 



2 In the theory of vortex dynam i cs for a two-dimensional ideal 
incompressible fluid (e.g. iLambl Il932t) . a conserved quantity 
playing the role of angular momentum is the angular impulse 
— J aj z r 2 dA. The torque required to move a small vortex 
patch of mass m and strength u> z and mass m radially at speed r 
is —muj z rr. This is a fraction — 2w z /f! of the torque required to 
move radially a particle of mass m in a circular Keplerian orbit 
of angular velocity fi. 
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situations involve the interchange of perturbed and ambient 
material. An equation that is ana logous to equ ation 14311 
holds in the case of convection (see ISchwarz schild 1953)- ^ n 
that case, the kinetic energy is proportional to the unper- 
turbed entropy gradient (analogue of the vortensity gradi- 
ent), gravity (analogue of the angular momentum gradient) 
and the square of the mixing length (analogue of the drift 
speed). 

Evidently this analysis may break down if the viscosity 
is too small, for then the anomaly AQ will not be small and 
further nonlinearity may intervene. In addition, if the vis- 
cosity if very small, this steady state may take so long to be 
achieved that the torque needs to be considered in a time- 
dependent sense. However, a possible way to rationalize the 
large value of t c obtained in the low-viscosity limit is as fol- 
lows. Suppose that the resonance drifts outwards in a disc in 
which £ is a decreasing function of r. By the time a steady 
state is reached, the ambient surface density, and therefore 
Tgt, ma y be very small. However, the torque may be sig- 
nificant because it derives from the trapped fluid. In this 
case t c will be very large, although it will not be described 
accurately by the approximations we have adopted. 

The second torque integral in equation 1411 . which de- 
rives from the open streamlines, is also a dimensionless quan- 
tity of order unity. In the limit of fast streaming (\v\ 2> 1) 
the dimensionless streamfunction and travel time may be 
approximated as x ~ — d8 + ^x 2 and A = —x/d. The second 
integral in equation (I4H then reduces to a standard Fresnel 
integral with the result that t c — * 1- The unsaturated GT79 
torque is therefore recovered in the limit of fast migration. 



5 ANALYSIS IN FOURIER SPACE 

As in Paper I, we apply a Fourier analysis in x and 0, writing 



«/ — oo „ „ 



dfe. (44) 



Equation (36) of Paper I then becomes 

G n 



— — ldk + n— + vk -77 

dt dk \l + k 2 



ka 



{Gn + l — G n -l) = (<5„,1 — 5n,~l) S(k). 



(45) 



The total tidal torque exerted in the corotation region (equa- 
tion 38 of Paper I) is, in dimensionless form, 



t c = Gi(0)-G-i(0). 



(46) 



Under the simplifying assumption b — 1 (Paper I) it is 
natural to rescale the wavenumber and time variable to 



k = v 1/3 k, U = 



so that 



d ~ d ,- \ 

— 2ivpk + n— + k )Gn 

dU F dk 



-pk (G n +l - G n -l) = (S n ,i - S n ,-i) 5(k), 



(47) 



(48) 



where 



p = -av 



-2/3 



-k 2 d In fi/dln r 



-m dfi/dr 



2/3 



(49) 



The dimensionless parameter p can be understood as the ra- 
tio of the characteristic time-scale of viscous diffusion across 
the librating region of the resonance to the characteristic 
time-scale of libration (raised to the power 2/3). In Paper I 
we found that p alone determines the degree of saturation 
of the resonance in the absence of migration. 

Linear theory applies when p< 1, although vp need not 
be small. Only modes n = ±1 are then significantly excited, 
and the solution in a steady state is 



G±i = ±H(±k)exp 



(50) 



giving rise to a dimensionless torque t c = 1, i.e. T c = Tgt. 
The corotation torque is therefore unaffected by planetary 
migration when p <C 1, although the form of the distur- 
bance is changed into a damped wave propagating radially 
downstream of the resonance. 



6 METHODS OF NUMERICAL SOLUTION 

We solve equation l|48^ and evaluate the corotation torque 
by two independent methods. The first method is identical 
to that described in Section 3.10 of Paper I, in which we 
obtain steady solutions by setting the time-derivative to zero 
and solving the (truncated) system of ordinary differential 
equations by a shooting method. 

The second method treats the equations as an initial- 
value problem, for which the initial condition correspond- 
ing to an unperturbed disc is G„(fc,0) = 0. The operator 
(d/dU) + n(d/dk) in equation (EHJ 

can be regarded as a 
Lagrangian derivative following a shearing flow in the semi- 
discrete (n, k) Fourier space. We therefore discretize the 
problem with respect to k by introducing a shearing lattice 
of Fourier wavevectors, 



k n ,j(U) — j 5k + nt* y 



(51) 



where j is an integer index and 5k is the lattice spacing. 
As the lattice shears, it periodically regains its form with 
a period of 5t, — 5k. The solution G is represented by its 
values G nj (i») at the moving lattice points, so that equation 
148H becomes 

{-ig %vpkn,j +k 2 ^ G n ,j 



~pk n j [G, 



+ij-(u/su) 



G 



n-l,j+(W<5<«)J 



= S nA 5(t, +jSU) - (5„_i<5(**- j6U), 
while the dimensionless torque is given by 



tc = G, 



-(t./«.) 



G_ 



l,(t« /St,)- 



(52) 



(53) 



In this representation, each Fourier amplitude starts at zero 
at t» — 0; the wavevector evolves according to the motion of 
the shearing lattice and the amplitude becomes non-zero 
either through the coupling terms (proportional to p) or 
through direct forcing (if n — ±1). Now t*/5t* is not gen- 
erally an integer, so some interpolation is required for the 
coupling terms; we use a simple linear interpolation between 
the bracketing integer values. 
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It is unnecessary to follow the wavevectors when they 
shear to large values of k because they become very strongly 
damped by viscosity. We apply a truncation such that \j\ < 
J, \n\ < N and periodically remap the wavenumbers, which 
is equivalent to resetting t* back to zero every time it reaches 
SU. 



7 NUMERICAL RESULTS 

The corotation torque in a steady state is shown in Fig. |5] 
For each value of the parameter p there is a certain degree of 
saturation of the torque at v = 0, as described in Paper I. As 
\v\ increases from zero the torque increases and can become 
significantly larger than the GT79 value (i.e. t c > 1) if the 
viscosity is small enough (i.e. if p is large enough). This 
enhanced torque derives from the librating region where a 
viscously equilibrated vortensity anomaly is established. As 
\v\ is increased further the torque is reduced because the 
librating region shrinks and disappears at \v\ = 1. Therefore 
the torque is maximized for \v\ fa 1/2. For |u| > 1 the torque 
slightly exceeds the GT79 value for any value of p. In this 
limit the torque derives from gas that streams through the 
resonance, and no saturation occurs. 

The general form of the numerically determined steady 
torque agrees fairly well with the approximate relation I4H . 
confirming our interpretation. It appears, however, that the 
convergence to the limiting form as the viscosity is reduced 
is rather slow. 

The time-dependence of the corotation torque, start- 
ing from an unperturbed disc, is shown for some repre- 
se ntative cases in F ig. El T he cas e v — is reminiscent 
of iBalmforth fc Korvcanskvl l)200lfl : the torque undergoes 
damped oscillations as the librating fluid stirs and mixes 
the vortensity, leading to saturation of the resonance. For 
small p, no oscillations occur and the degree of saturation 
is slight. For large p, many oscillations occur and the final 
steady torque is much reduced. 

The case v — 0.5 illustrates what happens when a sub- 
stantial asymmetrical librating region exists. The dynamics 
of that region no longer leads to a cancellation of the torque 
because of the systematic vortensity anomaly that builds up 
and is limited only by viscous diffusion. For large p, the time 
required to reach a steady state is considerable and the final 
torque is much enhanced. 

Finally, the case v = 2 corresponds to a situation of fast 
migration in which no librating region occurs. The distur- 
bance now takes the form of a kinematic wave downstream 
of the resonance. Care must be taken that the spatial do- 
main is large enough that the wave is viscously attenuated 
before it reaches the boundary. For large p some oscillations 
in the torque occur but the result differs only slightly from 
the GT79 value. Good agreement is found between the re- 
sults obtained with independent numerical methods. 



APPLICATION TO ECCENTRIC 
RESONANCES 



As discussed in Paper I and bv iGoldreich fc Saril {2003), 
the modification of the corotation torque has consequences 
for the eccentricity evolution of young planets orbiting in 



a protoplanetary disc. We consider the first-order eccentric 
corotation resonances associated with a planet of mass ra- 
tio q — Mp/M, executing an orbit of eccentricity e and 
semimajor axis a p within a Keplerian disc. Inner and outer 
resonances occur at radii 



\m±lj 



2/3 



(54) 



and the associated torques lead to eccentricity damping. The 
corresponding dimensionless parameters p and v for a mi- 
grating planet are 



p = 0.7006C^m 5/3 ega"^ /a ( — 



-2/3 



dr c /dt 



-4/3 



1.604C±m 2 eq V rQ, 



(55) 
(56) 



where the notation is as in Paper I (C„ being a correction 
factor that tends to unity for large m). Furthermore 



vp = 0.6552 



m \ r J 



1/3 



dr c /dt 



(57) 



where u„ = (3/2) (f/r) is the characteristic magnitude of 
the radial velocity in a Keplerian accretion disc. 

In type-II migration the factor (dr c /dt) /u„ is close to 
unity while the factor in square brackets is very small. Either 
v or p must be small, so the large enhancement of the torque 
seen in Fig. |2] does not occur. It is possible to have \v\ > 1 
if e is sufficiently small. This situation occurs because the 
librating region is narrow and the libration speed small, so 
the dimensionless drift speed is large. But then p <C 1, since 
the dimensionless diffusion time-scale is short, and so t c ~ 1 
whether the planet is migrating or not. Therefore planetary 
migration in the type-II regime is unimportant for eccen- 
tricity evolution, at least as far as the eccentric corotation 
resonances are concerned. 

In type-I migration the factor (dr c /dt)/tt„ can be much 
larger than unity and vp need not be small. Because of the 
dependence of v on m it is likely that the large enhancement 
of the torque is restricted to a few values of m at most. The 
issue is less important here because in the type-I regime 
the eccentricity damping is dominated by coorbital Lindblad 
resonances. 



9 TENTATIVE APPLICATION TO THE 
COORBITAL REGION 

9.1 Approximate flow in the coorbital region 

We now consider the possible application of these results 
to the more difficult problem of the coorbital region for the 
case of a planet of mass M p in a circular orbit of radius 
r p (t) undergoing radial migration. The perturbing potential 
in this case, including the indirect term, is 

-1/2 



$' = -GM p {r + r% - 2rr p cos ipY 
+GM p r~ 2 r cos ip, 

and the corotation radius is r c — (1 + q) 



(58) 
where q — 



M p /M* is the mass ratio. For ij < 1, and in the coorbital 
region where x — r — r c satisfies \x\ r c , the perturbing 
potential can be adequately approximated as 
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Figure 2. Steady dimensionless torque as a function of the dimensionless drift speed v for p = 0.5, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 (p = 0.5 is the 
flattest curve). The figure is symmetrical about v = 0. The approximate relation 1411 is shown as dotted lines for the same values of p. 



-GM p (r 2 s 2 + x 2 ) 



-1/2 



+ GM p r c cos ip, 



(59) 



where s = 2sin(^/2). The x-dependence in this expression is 
needed, of course, to localize the singularity of the potential 
at the location of the planet. In the case of the non-coorbital 
corotation resonance it was possible to treat the perturbing 
potential as being independent of x. 

We construct an approximation to the flow in the coor- 
bital region by taking the leading-order radial and azimuthal 
velocity perturbations to be 



1 d$' 



1 d& 
2Ql)x~' 



(60) 



These expressions derive from the linearized fluid dynamical 
equations in the corotation region, when terms proportional 
to x and those associated with pressure and viscosity are ne- 
glected. The expression for u' is identical to equation 1121 . 
but the equivalent for v' is also needed for consistency near 
the planet. Thus the leading-order velocity field in the co- 
moving frame is 



u = < 2grf2 sin tp 



2 . X 

S + —n 



+ 



-3/2 



-3/2 



— Q.x 
2 



(61) 



where r and Q are evaluated at r — r c , and we have spe- 
cialized to the case of a Keplerian disc. Fig. [1] shows some 
streamlines, based on this velocity field, for an inwardly mi- 
grating planet. There is a trapped region on the leading 
side of the planet, and there are open streamlines that pass 
the planet on its trailing side. This behaviour is a conse- 
quence of the higher (lower) magnitude of radial velocity 
that results on the trailing (leading) side of the planet when 
the inward migration velocity is combined with the veloc- 
ity (u',v'). A simi lar asymmetric behaviour was found by 
lArtvmowiczl (I2004T) . 

In the absence of migration, the stagnation points can 
be identified as the standard Lagrange points: 

L1.L2: ^t^ 1 '' 
L3: x = 0, 



if : 







L4, L5: 



0, 



if : 
If : 



(62) 



As the migration rate is increased from zero, a bifurcation 
occurs at |dr c /dt| « 1.45 qQ, at which what were L3 and L4 
(or L5, depending on the sign of dr c /dt) merge and disap- 
pear. This behaviour can be seen in Figs and |SJ where the 
separatrices are plotted. However, a librating region contin- 
ues to exist for any migration rate. This is an important 
difference with the case of the non-coorbital corotation res- 
onance. 
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Figure 3. Time-dependence of the dimensionless torque for various parameter values. The scaling of the time variable is related to the 
characteristic libration time (in the case v = 0, at least). Upper panel: v = and p = 0.4 (dashed), 2 (dotted) and 10 (solid). Middle 
panel: v = 0.5 and p = 0.4 (dashed), 2 (dotted) and 10 (solid). Lower panel: v = 2 and p = 0.1 (dashed), 0.5 (dotted) and 2.5 (solid). 



9.2 Stability of migration with a radial 
velocity-dependent torque 

If the coorbital torque depends on the migration rate o p , 
then dp satisfies a nonlinear equation of the form 

^d p = T c +T L , (63) 
aa p 

where J p = M p (GM*a p ) 1 ' 2 is the angular momentum of the 
planet in a circular orbit, and Tl is the 'Lindblad' torque 
that comes from the non-coorbital region of the disc and 
does not depend (or not sensitively, at least) on d p . If the 
migration is slow enough that the coorbital region is in a 
quasi-steady state, then T c can be replaced with its steady 
value T cs (d p ). 



In lMasset fc Papaloizoul (120031) it was argued that T c is 
directly proportional to d p , but evaluated at a retarded time 
owing to the effects of libration. Our analysis suggests that 
the dependence on d p is nonlinear in general, although in the 
coorbital case it may not resemble Fig. [5] in detail because 
of the persistence of the librating region as d p is increased. 
Equation 1631 may then have any number of solutions rep- 
resenting quasi-steady migration, and these solutions may 
be stable or unstable. 

Suppose that the planet migrates quasi-steadily accord- 
ing to one of the solutions of equation (1631 but is then per- 
turbed slightly, for example by applying and removing an ex- 
traneous torque. The corotation torque will take some time 
to adjust to the these changes. If it relaxes monotonically 
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Figure 4. Some streamlines in the coorbital region for the case of a planet with mass ratio q = 0.001 and inward migration rate 
dr c /dt = — 4<jf2. The figure is based on the approximation 1611 to the motion in the comoving frame with the planet located at the 
origin. The light streamlines are closed and trapped by the planet; the heavy streamlines are open. The trapped region is centred on the 
leading side of the planet and persists (but diminishes in size) at higher migration rates. The lowest streamline plotted is traced by the 
solid arrows, as it wraps in azimuth. This streamline circulates around the central star and then moves outwards, passing the planet, 
and continues outwards. 



towards the steady value T ca (d p ), then it can be shown that 
the migration is stable if 

p < (64) 
da p da p 

and unstable in the opposite case. [Formally this can be 
shown by combining equation Ki'M with a relaxation equa- 
tion such as dT c /di = (T cs — T c )/r.] If the quasi-steady 
solutions of equation 16311 are viewed as the points of inter- 
section of the graph r cs (fflp) with a straight line of gradi- 
ent dJp/da, then the solutions alternate in stability as the 
curve cuts the line alter nately from below and from above. 
While in the analysis of lMasset fc Papaloizoul ll200.?l) there 
is only one solution and it is either stable or unstable, here 
the possibility exists, in principle, of the planet finding an 
alternative migration rate that is stable. 



10 SUMMARY AND DISCUSSION 

We have investigated the effects of the radial migration of 
a planet through a disc on the corotation resonance and 
the associated torque. Our analysis deals mainly with non- 



coorbital corotati on resonances, which are crit ical to eccen- 
tricity evolution (|Ooldreich fc Tremaineiri98(t) . It was car- 
ried out using semi-analytical methods that account for the 
effects of gas pressure, turbulent viscosity and nonlinear sat- 
uration, as well the migration of the planet. Our approach 
does not account for the possible effects of shocks that cause 
changes in vortensity, which can be studied instead through 
nonlinear simulations. 

In the absence of radial motion, the corotational flow 
pattern consists of islands of closed streamlines, surrounded 
by a modified Keplerian flow that circulates around the star. 
The radial motion of the planet causes the resonance to drift 
through the disc and modifies the streamlines in the comov- 
ing frame, separating the closed regions and allowing radial 
flow through the coorbital region (see Fig.0. At faster mi- 
gration rates the closed regions are destroyed. The charac- 
teristic radial velocity for this transition is the radial width 
of the libration region divided by the libration time-scale. 

We find that the ratio of the resonant torque in a steady 
state to the value given by GT79 depends essentially on two 
dimensionless parameters. One of these (p) is proportional to 
the ratio of the characteristic time-scale of viscous diffusion 
across the librating region of the resonance to the charac- 
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Figure 5. Separatrices of streamlines in the coorbital region. The figure refers to the approximation 1611 to the motion in the comoving 
frame, for various values of the migration rate. The streamlines continuously fill the regions between separatrices. 



teristic time-scale of libration (raised to the power 2/3); in 
the absence of migration, this parameter alone determines 
the degree of saturation of the resonance. The second pa- 
rameter (v) is proportional to the ratio of the drift speed 
of the resonance to the characteristic radial velocity in the 
librating region; this determines the shape of the stream- 
lines. When the drift speed is comparable to the libration 
speed and the viscosity is small, the torque achieved in a 
steady state can be much larger than the unsaturated value 
in the absence of migration, but is still proportional to the 
large-scale vortensity gradient in the disc. 



The resulting corotation torque is a nonlinear function 
of drift speed, which peaks at the characteristic velocity (see 
Fig. The radial motion of a planet generally enhances 
the corotation torque by an amount that varies inversely 
with the turbulent viscosity, in the limit of small viscosity. 
This torque does not change sign under reversal of the direc- 
tion of radial motion. The corotation torque is proportional 
to the vortensity gradient, and the dependence on viscosity 
can be understood in terms of the vortensity redistribution 
in the corotation region. The trapped material attempts to 
preserve its original vortensity as it migrates radially into lo- 
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Figure 6. Continuation of Fig, l5l 



cations of different vortensity. As the contrast in vortensity 
between the trapped regions and surrounding background 
increases, the torque increases as well, but is limited by vis- 
cous diffusion between the two regions. If the viscosity is 
very small this steady state may take so long to be achieved 
that the torque needs to be considered in a time-dependent 
sense. A generally smaller contribution to the torque arises 
from material that passes across the coorbital region, be- 
tween the islands of closed streamlines. 

We suggest that the physics of the non-coorbital reso- 
nance can provide insight into the more complicated situa- 
tion for the coorbital region. There are similarities between 



the flow patterns in the two cases (Figs Q and 2J. In each 
case, there are regions of closed streamlines that move with 
the planet, as well as open streamlines that pass through 
the corotation region. We expect that the torque that re- 
sults from the closed streamlines in the coorbital case also 
involves the contrast between the background vortensity and 
the advected vortensity. The torque from the closed stream- 
lines would then be again inversely proportional to viscosity, 
in the small viscosity limit, provided that a large vortensity 
anomaly is able to accumulate. There are certainly differ- 
ences between the two cases. In the coorbital region, many 
different azimuthal components of the potential participate 
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and there is a singularity present due to the presence of the 
planet. Because of the singularity, there is always a region 
of trapped streamlines on one side of the planet for any mi- 
gration rate and this region will contribute to the coorbital 
torque. 

The analogy between the non-coorbital and coorbital 
cases also suggests that there may be, in some circum- 
stances, several possible migration rates in a given system 
because the corotation torque is a nonlinear function of the 
migration rate. Which one is achieved depends on the initial 
conditions and on considerations of stability, as discussed in 
Section 19.21 It is possible that fast migration does occur 
under some appropriate conditions. Mas set~fc Papaloizoul 
( 2003) emphasized the importance of a mass deficit in a par- 
tially cleared gap. Our analysis does not involve a gap, but 
considers the role of advected vortensity. The mass deficit 
would certainly contribute to the advected vortensity con- 
trast with the background, but may not be necessary. More 
accurate modelling of the coorbital region is required to de- 
termine its contribution to the migration rate. 
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APPENDIX A: ASYMPTOTIC SOLUTION OF 
THE FORCED ADVECTION-DIFFUSION 
EQUATION 



We consider the advection-diffusion equation 



+ u-VQ' -u\7 2 Q' 



(Al) 



dQ' 

at 

where u is a specified steady, two-dimensional velocity field 
satisfying V • « = 0, i/ > is a uniform diffusivity, and 
S is a specified steady source function. A streamline of u 
that does not include a stagnation point can either close 
on itself or be open to infinity. We assume that u has a 
region C of closed streamlines bounded by a separatrix and 
surrounded by a region O of open streamlines; the argument 
is easily generalized to allow for multiple closed regions. The 
boundary condition on open streamlines is that Q' — > 
far upstream. The interaction between the source and the 
response should be effectively confined to a limited region of 
space, either because S decays at large distance or because 
rapid phase variations lead to cancellation. 

Let x be a streamfunction such that u = Vx x e z , where 
e z is the unit vector normal to the plane of the flow, and let 
A be the travel time along a streamline (measured from its 
intersection with an arbitrary continuous curve). Then (A, x) 
are coordinates that cover C and O separately. In region C, 
A £ [0, A(x)) is a periodic variable on each streamline; in 
region O, A £ (—00,00). The element of area is AA = dAdx- 

A steady solution of equation (IAH satisfies 

^-u\7 2 Q' = S. (A2) 

We seek an asymptotic solution for large Reynolds number 
(small v). A naive expansion of the form 



Q' = Qo(A,x) + ^Qi(A,x) + - 

would imply, at O{u ), 
dQ' _ „ 



(A3) 



(A4) 



This generally fails in region C because it requires the source 
function to satisfy the solvability condition 







A(x) 



Sd\ 



(A5) 



on each streamline in C. Instead, the solution is of the form 
Q' = v-'Q'-iix) + Qo(A, X ) + fQi(A,x) + • ■ • . (A6) 
Equation 1A2II at O(v ) and 0(i/ 1 ) then gives 
dQ'o 



- v 2 o'_! = s, 



(A7) 



V 2 Qo = 0. 



(A8) 



d\ 

dQ'x 
d\ 

In region O the upstream boundary condition Q' — * as 
A — > — 00 requires that Q'-i vanish identically. The 'anoma- 
lous' response Q'_i is confined to the region of closed stream- 
lines on which there is a net forcing and the response builds 
up until limited by viscous diffusion. 
In region O we therefore have 



Q'o(A,x) 



5(A',x)dA'. 



(A9) 
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In region C the solvability conditions 

rHx) rHx) 

SdX, 



V 2 Q'_idA 



A(x) 



V 2 Qo dA = 



(A10) 
(All) 



apply on each streamline. Let A(x) denote the area enclosed 
by the (closed) streamline C(x) on which \ — constant. 
Then we have (after integration with respect to \) 



V 2 Q'_idA= / SdA, 

'A(x) Ja(x) 

The divergence theorem implies 



(A12) 



V 2 Q'-i d A 



Mx) 



c{x) 
dQ'_ x 



V<3'_! • n ds 

Vx ■ n ds 



c(x) 



V 2 X dA, 



-Mx) 



and so 
d X 



SdA 



(A13) 



(A14) 



Mx) 



'Mx) 



For continuity with region O in which Q'_i = 0, we have 
Q'-i(xs) = 0, where x — Xs is the separatrix. Therefore we 
have uniquely determined QLj in C and Q' in C It should 
be noted, however, that Q' as given by equation <A9t will 
be viscously attenuated at sufficiently large A. 
Now consider the integral 



Q'SdA 



(A15) 



over a region 1Z containing the part of the flow (including the 
whole of C) where the interaction is significant. This integral 
is directly related to the corotation torque. In a steady state 



J = 



]' (u ■ VQ' - vV 2 Q') dA 
/ (\q' 2 u-vQ'\7Q'\ -nds+ / u\\7Q'\ 2 dA. 



In terms of the asymptotic solution we find 
I = v' 1 !-! + Jo + ■ ■ ■ , 
with 

I-i= JjVQ'-x\ 2 AA, 

Io= I ^QoU-nds + 2 /(VQL0- VQ'odA 
Jan 1 Jc 

Although the advective flux in the first term in the expres- 
sion for Jo is converted into a viscous flux as the disturbance 
is attenuated, this term can be evaluated as if viscosity had 
no effect and the disturbance was advected to A — > +oo. The 
second term vanishes because 

JiVQ'^-VQ'o&A = - J^Q'_ x V 2 Q' dA 



(A16) 
(A17) 

(A18) 
(A19) 



V 2 Q dAd x 



= 0, (A20) 

where we use the boundary condition Q'_\ — on the sepa- 
ratrix and the solvability condition 1A111 . 
Thus the leading contribution 



J-i = 



f dQU 



|V X | 2 dA 



(A21) 



comes entirely from the closed region, while the next term 



Jo 





poo 


L 2 


/ SdX 

y - co 



dx 



(A22) 



comes entirely from the open region. 



